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Dynamical chaos in the problem of magnetic jet collimation 
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ABSTRACT 

We investigate dynamics of a jet collimated by magneto-torsional oscillations. The 
problem is reduced to an ordinary differential equation containing a singularity and 
depending on a parameter. We find a parameter range for which this system has stable 
periodic solutions and study bifurcations of these solutions. We use Poincare sections 
to demonstrate existence of domains of regular and chaotic motions. We investigate 
transition from periodic to chaotic solutions through a sequence of period doublings. 
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1 INTRODUCTION 

Many quasars and active galactic nuclei are connected 
with long thin collimated outbursts - jets. When observed 
with high angular resolution, these jets show structure 
with bright knots separated by relatively dark regions. A 
mechanism of collimation of such jets is still question- 
able. Magnetic collimation of jets was first c onsidered by 
iBisnovatvi-Kogan. Komberg fc FridmanI l|l969f ). In the pa- 
per of lBisnovatvi-KoganI ( 20071 ) magnetic collimation result- 
ing from torsional oscillations of a cylinder with elongated 
magnetic field and periodically distributed initial rotation 
around the cylinder axis was considered (Fig.[T]). The stabi- 
lizing azimuthal magnetic field is created here by torsional 
oscillations. An approximate simplified model was devel- 
oped, and an ordinary differential equation was derived de- 
scribing the process of dynamic stabilization. The interval 
of parameters, for jet stabilization to occur, was estimated 
qualitatively. 

The ordinary differential equation under consideration 
is a non-linear non-autonomous time-periodic second order 
equation with a singularity in the right hand side. The equa- 
tion contains one dimensionless parameter D, which summa- 
rizes the information about the magnetic field, amplitude 
and frequency of oscillations, radius of the jet, its spatial 
period along the jet axis, and sound speed in the jet matter. 
Here we investigate analytically and numerically the struc- 
ture of the phase space of this equation, which has a very 
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peculiar character and contains chaotic solutions as well as 
quasi-periodic and periodic regular solutions. 

The paper is organized as follows. In section 2 we de- 
scribe in more detail the mechanism of jet collimation by 
magneto-torsional oscillations and the dynamic confinement 
equation. In section 3 we discuss the main types of solutions 
of this equation. In section 4 we investigate analytically and 
numerically the dynamics of the system for large radii. In 
section 5 we construct Poincare sections for different values 
of parameter D. In section 6 we study periodic solutions 
undergoing a sequence of period doublings at values of pa- 
rameter D — Dn, n =1,2,3,. . . , when D increases. In section 
7 we discuss a possible mechanism of jet formation with a 
variable direction of angular velocity and estimate the value 
of D for astrophysical jets. 

We have found that stable periodic solutions disappear 
after an infinite cascade of period doublings. In the problem 
of period doubling the limiting constant q for values [D„ — 
Dn-\\/[Dn+i — Dn] is usually considered. For large n the 
constant q equals to Feige nbaum constant F = 4.6692... 
for dissipative systems, see iFeigenbaum (IIQSOT). and FH = 
8.721... for Hamiltonian systems, see lGreene et al.l (|l98ll ). In 
our work we obtain that q is approaching 8.72, in agreement 
with the expected behaviour for Hamiltonian systems. 



* E-mail: gkogan@iki.rssi.ru (GSBK); aneishta@iki.rssi.ru 
(AIN); seidov@bgu.ac.il (ZFS); tsupko@iki.rssi.ru (OYuT); 
krivosheev@iki.rssi.ru (YuMK) 



2 DYNAMIC CONFINEMENT OF JETS BY 
MAGNETOTORSIONAL OSCILLATIONS 

We consider jet stabilization by pure magneto- 
hydrodynamic mechanism associated with torsional 
oscillations. Jets are relativistically moving objects, and 
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Figure 1. Jet confinement by magneto-torsional oscillations 
(qualitative picture). 



the whole analysis of present paper takes place at the rest 
frame of the bulk motion of the jet. We suggest that matter 
in the jet is rotating, and different parts of the jet rotate 
in different directions (see Fig. [l|. Such distribution of 
rotational velocity produces azimuthal magnetic field, which 
prevents the disruption of the jet. The jet is represented by 
a periodic, or quasi-periodic structure along the axis, and 
its radius oscillates with time all along the axis. Space and 
time periods of oscillations depend on the conditions of jet 
formation: the length-scale, the amplitude of the rotational 
velocity, and the strength of the magnetic field. Time 
period of oscillations can be calculated in the framework of 
the dynamical model. The range of parameters, for which 
dynamical stabilization occurs, should also be inferrible 
from the model. Two-dimensional non-stationary MHD 
calculations are needed to solve the problem numerically. A 
ve ry simplified mode l of th is phenomenon was constructed 
bv lBisnovatvi-KoganI ([20071). This model allowed to confirm 
the possibility of such stabilization, to estimate the interval 
of parameters when it takes place and to establish the 
connection between time and space scales, magnetic field 
strength, and amplitude of rotational velocity. 

Let us consider a long cylinder with magnetic field di- 
rected along its axis. This cylinder will expand without 
limit under the action of pressure and magnetic forces. It 
is possible, however, that a limiting value of the cylinder ra- 
dius could be reached in a dynamic state, when the whole 
cylinder undergoes magneto-torsional oscillations. Such os- 
cillations produce a toroidal field (magnetic field lines are 
frozen into matter), which prevents radial expansion. There 
is therefore a competition between the induced toroidal field, 
compressing the cylinder in the radial direction, and the 
gas pressure, together with the field along the cylinder axis 
(poloidal), tending to increase its radius. During magneto- 
torsional oscillations there are phases when either compres- 
sion or expansion force prevails, and, depending on the in- 
put parameters, there are three possible kinds of behaviour 
of such a cylinder that has negligible self-gravity. 

(1) The oscillation amplitude is low, so the cylinder suf- 
fers unlimited expansion (no confinement). 

(2) The oscillation amplitude is high, so the pinch action 
of the toroidal field destroys the cylinder and leads to the 
formation of separated blobs. 

(3) The oscillation amplitude is moderate, so the cylin- 
der, in absence of any damping, survives for an unlimited 
time, and its parameters (radius, density, magnetic field, 
etc.) change periodically, or quasi-periodically, or chaotically 
in time. 

These phenomena can be described in general 
by the system of ax ially symmetric MHD equations 
l|Bisnovatvi-Koganll2007l ) . The system of equations was sim- 



plified in the paper of iBisnovatvi-KoganI (|2007f ) to investi- 
gate the most important property of dynamical competition 
between different forces in order to check for the possibility 
of dynamical confinement. A profiling procedure was used 
for this purpose. In particular, gravity in the direction of the 
cylinder axis was neglected, approximate uniform density 
along the radius was assumed, and adiabatic case with poly- 
tropic equation of state was considered. Approximate system 
allowed to investigate linear oscillations (around the equilib- 
rium state in presence of gravity) of infinite, self-gravitating 
cylinder with uniform magnetic field and rotation along its 
axis. 

For further simplification gravity was neglected, be- 
cause in a relativistic jet the self-gravitating force is ex- 
pected to be much lower than the magnetic and pressure 
forces. In this approximation we have the following feature 
of problem. Without gravity the equilibrium static state of 
the cylinder does not exist. But it turns out that neverthe- 
less bounded solutions exist, where cylinder radius oscillates 
and remains finite (dynamic confinement). Confinement by 
magneto-torsional oscillations is therefore physically realiz- 
able. Similar situation is observed for instance in the sim- 
ulations of electron-positron discharges in the polar cap of 
pulsar - the pair cascade settles down to a stable state whe re 
it has a limiting cycle type of behaviour l|Timokhinll2010^ . 

After considerable simplifications, details of which may 
be found in the paper of iBisnovatyi-K ogan (2007), the equa- 
tion, describing the magneto-torsional oscillations of a long 
cylinder, takes the following form: 



dr2 



1 - D sin^ T 

y 



(1) 



This equation describes approximately time dependence of 
the outer radius of the cylinder R{i) in the symmetry plane, 
where the rotational velocity remains zero. The cylinder has 
magnetic field Bz it) , isothermal equation of state of matter 
P = Kp, maximal amplitude of the angular velocity of os- 
cillations f2o, and angular frequency of oscillations uj. The 
oscillating cylinder has a periodic structure along the z axis, 
with space period zq. The nodes with zero amplitude of os- 
cillations are situated ai z = ±n^, and the oscillations with 
maximal angular amplitude fio are situated at z = -^zbn-^, 
n = 0, 1, 2.... Axial motion of the matter in the cylinder is 
neglected, and density p{t) is supposed to be uniform along 
the radius. Under these conditions the equations for conser- 
vation of mass and magnetic flux (we assume infinite elec- 
trical conductivity) determine the constant values Cm and 
Cf. 



p R — Cm = PoRo; Bz R — Cb = BzfiRo^ 



(2) 



where po, Rq, and B^^ are some characteristic values. The 
dimensionless variables and the parameter _D in ([TJ are de- 
fined as 



T = ujt, y = -^, 

with Ro ~ 



(3) 



D = 



_i / Cfc^o y 

KCm V ZqOJ J 



U) 2lYKCm \ ZqOO 

The frequency of oscillations u may be represented as 



U! = ankVA 



(4) 
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Figure 2. Time dependence of non-dimensional radius y (upper 
curve), and non-dimensional velocity z = y' (lower curve), for 
D = 2.0. 



Figure 3. Time dependence of non-dimensional radius y (upper 
curve), and non-dimensional velocity z = y' (lower curve), for 
D = 2.4. 



where k is the wave number, k — 2-k/zo, and Va is the 
Alfven velocity, Va = Bzfi/\/4:npo; a„ < 1 is a coefficient 
determining the frequency of non-linear Alfven oscillations, 
which are identical to magneto-torsional oscillations under 
investigation. In our problem there are no static solutions. 
In a balanced non-compressible cylinder the frequency of 
magneto-torsional oscillations is defined by Q with a„ ~ 1. 



3 GENERAL TYPES OF SOLUTIONS OF 
DYNAMIC CONFINEMENT EQUATION 

Equation ffl wa s solv ed numerically in the paper of 
iBisnovatvi-KoganI (|2007l ) for different values of parameter 
D, and fixed initial conditions y{0) — l,y'{0) — 0. The fol- 
lowing three regimes of the behaviour were found in these 
calculations. 

(1) D < 2.1 The oscillation amplitude is low, so the 
cylinder should suffer unlimited expansion (no confinement). 
This regime corresponds to small D. For example at D = 2 
there is no confinement, and radius grows to infinity after 
several low-amplitude oscillations. Fig. [21 

(2) D > 2.28 The oscillation amplitude is too high, so 
the pinch action of the toroidal field destroys the cylinder, 
and le ads to formation of sepa rated blobs. The calculations 
show l|Bisnovatvi-Koganll2007^ . that at D = 2.28 and larger, 
the radius finally goes to zero with time, but with different 
behaviour, depending on D. At D between 2.28 and 2.9 time 
dependence of the radius y may be very complicated, consist- 
ing of low-amplitude and large-amplitude oscillations, which 
finally decay to zero. The time at which radius becomes zero 
depends on D in a rather peculiar way, and it may happen 
at T < 100, like at D = 2.4 (Fig.[3ll, 2.6, or at r ~ lO'' like 
at D = 2.5 (in the last case the radius passes through very 
large values and then returns back to zero). For D — 3 and 
larger the solution is very simple: the radius goes to zero at 
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Figure 4. Time dependence of non-dimensional radius y (upper 
curve), and non-dimensional velocity z = y' (lower curve), for 
D = 3.0. 



T < 2.5, before the return of the right hand side of ([T]) to a 
positive value. Fig. (4] 

(3) For the intermediate interval 2.1 < D < 2.28 the 
radius is not growing to infinity, but is oscillating around 
some average value. The cylinder survives for a unlimited 
time, and its parameters (radius, density, magnetic field, 
etc.) change periodically or quasi-periodically in time. Figs 

ElEl 

The above results from the paper of iBisnovatvi-Koeanl 

(J2007t) have been obtained for solutions only with initial 
radius yo = 1- The solutions at moderate D are not pure 
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Figure 5. Time dependence of non-dimensional radius y (upper 
curve), and non-dimensional velocity z = y' (lower curve), for 
D = 2.1. 




Figure 6. Time dependence of non-dimensional radius y (upper 
curve), and non-dimensional velocity z = y' (lower curve), for 
D = 2.25. 



periodic, but it is not clear from the figures, whether are 
they regular or chaotic. 

In the current paper we investigate the solutions of ([T|) 
for different values of D and different initial radii j/o, mainly 
for moderate D's. Solutions of principal interest are those 
that do not go neither to nor to infinity. Such solutions 
correspond to stabilized jets. We use averaging method and 
Poincare sections to investigate the properties of these solu- 
tions. 




Figure 7. Phase curves of the system with Hamiltonian (|9j for 
D= 1.5, and W = -0.1,-0.05,0,0.05,0.1. 




Figure 8. Phase curves of the system with Hamiltonian ((9)l for 
D = 2.5, and H = -0.1, -0.05, 0, 0.05, 0.1. 



4 DYNAMICS FOR LARGE RADII 

Let us denote z — y' = dy/dr. Equation ([!]) is equivalent to 
the system 



y = z, z = 



D sin^ r 

y 



(5) 



We consider a motion for large values oi y. In this limit let 
us introduce a small number e and denote q — ey. We will 
assume that y ~ 1/e, and therefore g ~ 1. The equation 
system for q, z is 



q = ez, z = e- 



D sin T 



This is a Hamiltonian system with Hamiltonian 



eH 



(1 — Dsin r)lng 



and equations of motion 



g' = edH/dz, z — —edH/dq. 



(6) 



(7) 



(8) 



For small e the variables g, z are changing slowly 
with respect to changing of r. The averaging method 
(jBogolvubov fc Mitropolskiil Il961f) prescribes to average 
Hamiltonian (O over the fast variable r for an approximate 
description of the behaviour of the slow variables g, z. We 
get the averaged Hamiltonian 
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Figure 9. Poincare section near tlie stable fixed point, z = y' = 
dy/dr, D = 2.001. One periodic and five quasi-periodic solutions 
are shown. This figure is constructed numerically by solutions of 
equations Jsjl. Difference between curves in this figure and phase 
curves from approximate formula (I12I I is about thickness of line 
on figure, therefore we do not show the latter curves here. 




Figure 11. j/-coordinate of the stable fixed point as a function of 
D. In this figure the separated points are obtained by numerical 
calculations, the curve is plotted analytically with formula (llll l. 



0,02 
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Figure 10. Poincare section for D = 2.001. This figure is con- 
structed numerically by solutions of equations Jsjl. Difference be- 
tween curves in this figure and phase curves from approximate 
formula l|12p is about thickness of line on figure, therefore we do 
not show the latter curves here. 
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Figure 12. Phase trajectory corresponding to periodic solution 
for D = 2.001; initial data y{0) = 15.8034860625, 2(0) = 0. Phase 
trajectory is constructed by numerical solution of equations JSjl 
at given initial values and putting all dots of solution on the plane 



eU 



(i-f)i... 



(9) 



y 



Phase curves of the system with Hamiltonian (|9} for D < 2, 
and for D > 2 are shown in Figs [7] and [8l respectively. 

li D < 2, then for Hamiltonian Q q{t) and y{t) — >• 
±oo,z(t) — >■ ±00 as t — !> ±cx) for all solutions of the av- 
eraged system. This situation corresponds to the unlimited 
expansion of the jet (no collimation). If D > 2, then q{t) and 
y{t) — > 0,z{t) — > ±oo as i — >■ =F''"», where r, is some finite 
moment of time that depends on the chosen solution. This 
corresponds to strong initial amplitude of pulsations leading 
to the fragmentation of the jet into separate clumps. 

Solutions of the averaged system describe approxi- 
mately the solutions of the exact system for large values 
of y provided that {D — 2) is not small. 

It turns out that for small (D — 2) > the exact system 
has a stable periodic solution, which is described approxi- 
mately by the formulas 

1 1 



2^/D~~2 2 



VD - 2 cos 2r, z = VD-2sm 2r. (10) 



(see Appendix). In the Poincare section {r = Omod tv} this 
solution is depicted as a fixed point of the Poincare return 
map. This fixed point is located on the axis z — 0, its in- 
coordinate is given approximately by the formula 



2/D^^ 2 



vir^. 



(11) 



This fixed point is surrounded by a family of invariant curves 
of the Poincare return map. The approximate equation for 
these invariant curves is 



1,2 1 

2^ 16t/2 



4- 



2) In 2/ = const 



(12) 



(see Appendix). The fixed point and surrounding it invariant 
curves are shown in Fig. [9] for D = 2.001. 

The Poincare sections are constructed numerically in 
this paper. For the same value of D, we solve equations ((SJ 
for different initial j/o and zq. During the numerical integra- 
tion of equations, we mark moments t — mod n. For each 
integration, we put the points on the plane {y, z) at the mo- 
ments r = mod tt. The regular oscillations are represented 
by closed lines (invariant curves) on the Poincare maps: 
closed curves correspond to regular quasi-periodic solutions, 
points in the centers of families of such curves correspond 
to stable regular periodic solutions. Chaotic behaviour fills 
regions of finite area with dots on the Poincare maps. 

For this value of D = 2.001 the family of invariant 
curves covers rather a big domain in the Poincare section, 
see Fig. 1101 For initial conditions on the invariant curves the 
solution of the equation ([l]) is represented by quasi-periodic 
functions of r. The values of j/-coordinate of the fixed point, 
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D = 2.15 








Figure 15. Invariant curves on Poincare section for D = 2.25. 
The stable fixed point in the centre has coordinates y = 0.887, 
^ = 0. 



Figure 13. Poincare section for D = 2.15. Central point corre- 
sponds to a stable periodic solution of period tt. 



D = 2.15 




Figure 14. Zoom of the part of Fig. 1131 islands on Poincare 
section for D = 2.15. 



obtained theoretically (see (|11|) ') and numerically, are shown 
in Fig. [Til 

One can see that the obtained asymptotic expression 
for the coordinates (radii) of the stable fixed points works 
very well even for D not very close to 2 (and therefore when 
y is not very big). The projection onto the plane y,z (in 
another words, phase trajectory) of the periodic solution 
corresponding to the fixed point in Fig. [9] is shown in Fig. 



m 



5 POINCARE SECTIONS FOR DIFFERENT 
VALUES OF PARAMETER D 

In order to find bounded regular solutions we have con- 
structed a Poincare section for the system ((S)) for r = 
mod TT for different values of parameter D. We use such 
construction of the Poincare section because the right-hand 
side of the equatio n ([T]) is a periodic function of t with period 
TT. Recall that in iBisnovatvi-KoganI ((20071) such solutions 
with yo = 1 were found to exist only for ~ 2.1 < D <~ 2.28. 
We will concentrate here on the behaviour of the invariant 
curves on the Poincare section around the stable fixed point 
of the Poincare return map, and on the transition to chaotic 
behaviour in this region with increasing of parameter D. 

In Figs 1131 1141 the Poincare section for D = 2.15 is 
presented. The closed invariant curves cross axis z = in 
the interval 0.8 ~< yo <~ 1.80. Outside this interval solu- 
tions of equations ((5]) show non-regular chaotic behaviour, 
with unbounded trajectories. Central point in Fig. [13] cor- 
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Figure 16. Appearance of chaos on Poincare section for D 
2.25. 



responds to the periodic solution with period tt. One can 
see a chain of seven stability islands, centres of which deter- 
mine periodic solution with period ln\ zoom of this chain is 
presented in Fig. 1141 In Fig. [13] one can also see a chain of 
fifteen stability islands, centres of these islands correspond 
to periodic solution with period IStt. 

Invariant curves around the periodic solution for D — 
2.25 are plotted in Fig. 1151 In Fig. [16] the appearance of a 
chaotic solution is shown. The last shown invariant curve 
corresponds to the solution with initial values j/(0) — 0.757, 
z(0) = 0, while the solution with y{0) = 0.756, z(0) = is 
already chaotic. In Fig. [T7] a general picture of regular and 
chaotic trajectories is represented. 

Invariant curves and chaotic trajectories on the 
Poincare section are plotted for D — 2.30 in Fig. 1181 and 
for D = 2.35 in Fig. \T9\ 

The case D = 2.38 is represented in Fig. [20] Layers with 
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Figure 17. Invariant curves and chaotic trajectories on Poincare 
section for D = 2.25. 
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D = 2.30 




Figure 18. Invariant curves and chaotic trajectories on Poincare 
section for D = 2.30. Central point corresponds to a stable peri- 
odic solution of period it. 



D = 2.35 




ri,m 0.6S 0.70 0.72 qj4 y oje 0.78 aso o.s2 

Figure 19. Invariant curves and chaotic trajectories on Poincare 
section for D = 2.35. Central point corresponds to a stable peri- 
odic solution of period tt. 




Figure 21. Time dependence of non-dimensional radius y (upper 
curve) and non-dimensional velocity z = y' (lower curve) for pe- 
riodic solution of period tt, see central point in Fig. 1201 D = 2.38. 




Figure 22. Phase trajectory in the plane y, z for periodic solution 
of period tt, see central point in Fig. 1201 D = 2.38. 



both regular and chaotic trajectories are situated around pe- 
riodic solutions. In Figs l21|[22l tinie evolution j/(r), z{t) and 
phase trajectory in the plane y, z are shown for a periodic 
solution with period tt. 






ilS^^ivi':- 



D = 2.38 



i^y- 




Figure 20. Invariant curves and chaotic trajectories on Poincare 
section for D = 2.38. Central point corresponds to a stable peri- 
odic solution of period tt. 



6 TRANSITION TO CHAOS VIA CASCADE 
OF PERIOD DOUBLINGS 

With increasing D the stable periodic solution considered in 
Sections 4 and 5 loses its stability via period doubling bifur- 
cation. The first period doubling occurs at D\ — 2.38101. At 
this value Di the solution with period n becomes unstable, 
and a stable periodic solution of period 2n appears. 



D = 2.33100 




0,676 0,677 y 0,678 0,67S 

Figure 23. Poincare section for D = 2.38100. 
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D = 2.38101 
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Figure 24. Poincare section for D = 2.38101. This value of D 
corresponds to first period doubling. The left central point maps 
to the right central point after time vr, then visa versa. Compare 
with Fig. 1231 where Poincare section just before this doubling is 
shown. 



z 

0,0Q6 



D = 2.3875 
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Figure 25. Poincare section after first period-doubling; D = 
2.3875. Two points in centres of concentric curves correspond to 
the periodic solution of the period 2tt. 




Figure 26. Time dependence of non-dimensional radius y (upper 
curve) and non-dimensional velocity z = y' (lower curve) for the 
periodic solution of period 2tt, see points in Fig. 1251 D = 2.3875. 
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Figure 27. Phase trajectory in the plane y, z for periodic solution 
of period 2n, see points in Fig. [25] D = 2.3875. 



D = 2.39 





Figure 28. Poincare section for D = 2.39. Two points in the 
centres of concentric curves correspond to a periodic solution of 
period 2it. 



The subsequent points of period doubling Di, D2, ... 
can be found by constructing Poincare sections for different 
values of D. In Fig. [23] the central part of the Poincare sec- 
tion for value of D before the first bifurcation is shown. In 
Fig. [24] the central part of the Poincare section for value D 
soon after doubling is shown. Periodic solution with period 



D = 2.39 




tl-'.i 



0,708 0,710 0,712 0,714 y 0,716 0,718 0,720 

Figure 29. Zoom of left part of right island in Fig. 1281 The island 
at the Poincare section for D = 2.39. The point in the centre of 
the island corresponds to the periodic solution of the period 127r. 
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■K is now unstable, but there is a stable periodic solution of 
period 2-k. On Poincare section this solution is represented 
by two points surrounded by closed curves (Fig. [24]). Under 
one iteration of Poincare return map the left point and sur- 
rounding it curves are mapped, respectively, onto right point 
and surrounding it curves, and vice versa. With increasing 
D the distance between the points of periodic solution on 
Poincare section increases. 

In Fig |25l the Poincare section for a slightly greater value 
of D, namely D = 2.3875 is shown. In Figs [26l [27] time 
evolution j/(r), z{t) and phase trajectory in y,z plane for 
the periodic solution with period 27r are shown. 

At D2 = 2.410172 the periodic solution with period 
27r becomes unstable, and a stable periodic solution of pe- 
riod 47r appears. The following doubling points found by the 
Poincare section construction are: 



D^ 


~ 2.38101, 


TT — >• 2tT, 


D2 


~ 2.410172, 


2-K — >• 47r, 


Di 


~ 2.4137115, 


4n — ^ Stt, 


Di 


~ 2.4141168, 


Stt -^ 16tt, 


Ds 


~ 2.4141633, 


IBtt -> 327r 



The values Di, D2, D3, D4, have also been 
obtained with the help of of AUTO-07P: Con- 
tinuation and Bifurcation Software for Or dinar y 

Differential Equations (JDoed el & O ldemanl (|2009l l. 
http://cmvl.cs.co ncordia.ca/auto/ ). 

According to JGreene et al.l (JIQSIJ ) the limiting constant 



D„-D„ 



D 



n + l 



D„ 



(13) 



should for large n approach the constant FH — 8.721 for 
Hamiltonian systems. We obtain the following values: 



gi23 



?234 



D2-D1 






D2 
D2 



?345 = 



D4 



D3 
D3 



D5-D4 



8.24, 
8.73, 
8.72. 



We conclude that this behaviour agrees with the expected 
one for Hamiltonian systems, despite of the additional sym- 
metry z — >■ ~z, t — > —t. 

We can find approximately the value of parameter Doo 
for which stable periodic solutions disappear after an infinite 
cascade of period doublings using the value q ~ 8.721: 



Doc = Di + {D5 - D4) + (De, - D5) + {D7 - De) + ... 

= Di + (D..^Di) + £^-£± + £^^ + ... 

q q^ 

= D4 + •^'^•^'^ ~ 2.4141693. (14) 

1- 1/g 

In Table 1 we summarize the figures of the present ar- 
ticle. 



Table 1. Values of parameter D and numbers of corresponding 
figures. 



Values of parameter D 


Numbers of figures 


1.5 


m 


2.0 


m 


2.001 


MMWM 


2.1 


m 


2.15 


mm 


2.25 


|6] [15] [16] [17] 


2.30 


m 


2.35 


m 


2.38 


[20] [21] [22] 


2.38100 


m 


2.38101 


m 


2.3875 


[25] [26] [27] 


2.39 


mm 


2.4 


m 


2.5 


m 


3.0 


m 



7 DISCUSSION 

In this section we would like to address the issue of a jet 
angular velocity oscillations generation. Jets are associated 
with accretion and consist of ejected matter from accretion 
disc. Jet matter inherits angular momentum from the accre- 
tion disc matter. 

The jet matter may have a periodically (or quasi- 
periodically) varying sign of the angular momentum of the 
accretion disc matter around a supermassive black hole. In 
a dense stellar cluster around the SMBH, the accretion disc 
can change its direction of rotation due to different sign of 
the angular momentum of matter of the disrupted accret- 
ing stars. In this situation the magneto-torsional oscillations 
should be inevitably generated in the outflowing jets. 

In the case, when stellar cluster is rotating, the dis- 
rupted stars will preserve the excess of the angular momen- 
tum, and the jet may rotate, and there should be oscilla- 
tions of the angular velocity around the regular rotation. 
The presence of these oscillations will lead to effects sim- 
ilar to ones described in our model. But the presence of 
regular rotation makes the physical picture of phenomenon 
more complicated, because of the centrifugal force created 
by regular rotation and leading to jet expansion in radial 
direction. 

In the framework of our model we have demonstrated 
that for a narrow range of the parameter D the cylinder 
radius remains finite for a long time. This can potentially 
lead to a long-live jets. Due to complicated nature of the 
problem, we can not estimate the number of jets coUimated 
by this mechanism. In the case of D < 2 the jet collimation 
does not occur. We can tell that if D > 2 then jet exists in a 
continuous form or in the form of blobs, both configurations 
corresponding to observations. 

We can estimate value of D 



D 



1 



\ ZoijJ J 



2ttKC„ 
in the following way: 



no 



Rq 



Zo : 



7?0, Cn 



' poRo, 



(15) 



Cb ~ B.fiRi, (16) 
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K: 



c 



(17) 



Here c is the speed of light, _Ro, po and B^^ are some char- 
acteristic values of length, density and magnetic field corre- 
spondingly, A is the non-dimensional constant. We obtain: 



D = 



2n^al 



(18) 



For estimation of A we can write: P = nkT, n = -22- where 
P is the pressure, n is the concentration, T is the temper- 
ature, nip is the proton mass, k is the Boltzmann constant. 
Comparing with P — Kpo we obtain: 



A^ 



kT 



(19) 



For estimation in a nonlinear regim e of os cillations we can 
put al = 0.1 (see iBisnovatvi-Kogaiil (|2007l )'). 

We see that for a sufficiently large proton loading {A > 
1) we can have D > 2. 
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APPENDIX A: APPROXIMATE FORMULAS 
FOR PERIODIC SOLUTIONS AND INVARIANT 
CURVES 

We will use a standard approach for the averaging method. 
We will make a canonical time-periodic transformation of 
variables such that the Hamiltonian for the new variables 
will not depend on time in the principal approximation. Dis- 
carding small time-depending terms in the Hamiltonian we 
will get a Hamiltonian for some autonomous system. We 
will find an equilibrium position of this system. The exact 
system in the new variables has a periodic solution close 
to this equilibrium position. Then we will find formulas for 
this periodic solution in the old variables making use of the 
formulas for the transformation of variables. 

We will construct the required transformation of vari- 
ables as a composition of several transformations of vari- 
ables. In the system with Hamiltonian (O we make a canon- 
ical transformation of variables (g, z) i— >■ {q, z) with a gener- 
ating function of the form 



qz + eSi{q,z,T). 



(Al) 



The old variables (g, z) and the new variables (g, z) are re- 
lated via formulas 



z = z + e- 



dSi 

dq 



= q + e- 



dSi 
dz 



(A2) 



The dynamics of the new variables is described by the Hamil- 
tonian 



rr TT dSl 

eH = en+e-- — —e 
or 



1 f_^ dSi 



;-?(- 



cos2r) ) lng+ — — 



(A3) 



Let us choose Si such that the terms of order e in this Hamil- 
tonian will be independent on r: 



— - In g cos 2r + — — = 
2 OT 



One can 


take 






In g sin 2t. 


Then 




z = z + 


e — sin2T, q — q. 
4g 


The new 


Hamiltonian is 


eH = 


e 


i(^.--f.|sin2.j +\(D-2)lr.q 


= 


e 


1 2 Dz 2 D'^ 2 


+ 


\{D~2)\nq 





(A4) 
(A5) 
(A6) 



(A7) 



Now we make a canonical transformation of variables 
(g, z) I— >■ (g, z) in the system with Hamiltonian eS with the 
generating function of the form 



qz + e S2{q,z,T). 



(A8) 



The old variables {q, z) and the new variables (q, z) are re- 
lated via formulas 



z + e 



dq 



q = q + e 



>dS2 



(A9) 



The dynamics of the new variables is described by the Hamil- 
tonian 



eH 



eH + £■ 



■dS2 
dr 



1 /-^ 2dS2 



D f.^ 2dS2 , . „ 
- £-— z + £ -—- I sm 2r 
4(j \ dq 



32-,sm^2. + -(D-2)ln, + .— 



Let us choose S2 in such a form that the terms of order e 
in this Hamiltonian will be independent on r: 



Dz . ^ ^dS2 „ 
-— sm2r + — — = 0. 
Aq or 

One can take 



02 = — — cos 2t. 
Sq 

Then 



2DZ _ 2-D 

z — z — £ ^— ■ COS 2r, g = g + e ^ cos 2r. 

The new Hamiltonian is 



(AlO) 



(All) 



(A12) 



£H = £ 



If- 2DZ 2 sD^i 

— (z — £ — T- COS 2r) — e 

2^ 8g2 ' -iOn-i 



32g3 



COS 2r sin 2r 



£^D'^ 1-cos2t 1,^ 
+ 32i^^^ + 2(^-2) In g 



(A13) 



In this relation one should express g through g. Let us 
suppose that (D — 2) ~ e^ (only for such values of D the 
periodic solution will exist). Then we will just replace g with 
g; this will lead to the error 0(£^) in the Hamiltonian. Thus 
we have 
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that this system has an equilibrium position. We should cal- 
culate partial derivatives of the Hamiltonian and find a point 
where they vanish. Thus for the coordinates of the equilib- 
rium position we have 



5 = 0, 



£^D 



32g 
Thus we have 



.3+2i(^ 



2) =0. 



16(D-2)' 



(A16) 



(A17) 



In the numerator of this formula we can put D = 2 as 
this will lead to an error 0(£^) in g. Thus for an approximate 
position g, of the equilibrium we have 

'' - i7#=1- (^18) 

Plugging into previous formulas for transformation of vari- 
ables the coordinate of the equilibrium {q,z) = (g, ,0) we 
find a periodic solution of the original system. Approximate 
formulas for this solution are 



2VO-2 



2 D 

- £ - — COS It ! 



> 2 COS 2r 



8g. 



2VO-2 



-\/D-2cos2r, 



2^/D-2 

D 2 sin 2r r- 

z « £ sin2r = e ;^^^= = V-D — 2sin2r. 

4g. 4,£/2^/Tr^ 

Returning to y = g/e we get 

1 1 



V 



2^/5~^ 2 



VD - 2 COS 2r, z = V-D-2sin2r. (A19) 



The function H is an approximate first integral of motion: 
H « const. From this relation we obtain an approximate 
expression for the invariant curves of the Poincare return 
map: 



1 2 1 1 

2^ ^ 1%2 + 2^^"^^^°^"°°°''*' 



(A20) 



£H = 



£ 


[1,2 2D~Z^ 
2'' ^ 8g2 


\{D~2)\uq 


H 



cos 2r 



+ 0(£* 



£^D^ l-cos2r 
32g2 2 



(AM) 



Now we make a canonical, 0(e^)-close to the identical trans- 
formation of variables (g, z) 1— > (g, z) in the system with 
Hamiltonian £H such that the terms of order e"^ in the new 
Hamiltonian will be independent on r. The new Hamilto- 
nian in the principal approximation is just the average of 
the old Hamiltonian over r. Thus the new Hamiltonian is 



£H = £ 



1 p D 1 

-i' + ^-^ + -(-D-2)lng 

2 64g2 2 ^ ' ^ 



+ 0{£'). 



(A15) 



Let us neglect the term 0{£*') in this Hamiltonian. We obtain 
an autonomous Hamiltonian system. Let us demonstrate 



